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Abstract 

We performed two-dimensional radiation hydrodynamic simulations of supercritical accretion flows 
around neutron stars (NSs). In contrast with the accretion flows onto black holes (BHs), we find that 
the shell-shaped high-density regions form around the NSs, since the radiation force is enhanced in the 
innermost regions. The enhanced radiation force drives strong outflows above and below the disk. The 
mass-accretion rate onto the NS exceeds the critical rate, Le/c 2 , with Le being the Eddington luminosity. 
However it is about 20 — 30% of that onto the BH, under the condition that we employ the same mass-input 
rate, Mi nput , which is mass injected from the outer disk boundary per unit time. The mass-outflow rate 
is a few-times larger in flows around NSs than in flows around BHs. The supercritical NS accretion flows 
mainly release the accretion energy as the kinetic energy of the outflows, though the disk luminosity is 
predominant over the kinetic energy output rate in the BH accretion flows. The resulting velocity and 
mass-outflow rate of the outflows are 0.2 — 0.3c and 150 — 700Le/c 2 , respectively, for the mass-input rate 
of 3 x 10 2 <Af in p U t/(i E /c 2 )<3 x 10 3 . This implies that the SS433 jets can be roughly explained by the 
supercritical accretion onto a NS. However, the collimation angle of the outflows in our simulations (~ 20°) 
is larger than that of the SS433 jets (a few degrees). 

Key words: accretion, accretion disks — ISM: jets and outflows — hydrodynamics — stars: neutron 
— stars: individual (SS433) 



1. INTRODUCTION 

Luminous compact objects are thought to be powered 
by accretion flows. It is widely believed that accretion 
processes can be basically described by the standard disk 
model, proposed by Shakura & Sunyaev (1973). However, 
this model breaks down if the mass-accretion rate is com- 
parable to, or exceeds, the critical rate, Le/c 2 , where Le 
is the Eddington luminosity. Then, advective cooling is 
predominant over radiative cooling. 

Abramowicz et al. (1988) established the slim-disk 
model, which is one-dimensional model that includes 
the advective cooling. Two-dimensional radiation- 
hydrodynamic (RHD) simulations of supercritical accre- 
tion flows around black holes (BHs) were attempted by 
Eggum, Coroniti, & Katz (1988) and Okuda (2002). By 
performing long-term two-dimensional RHD simulations, 
we have for the first time confirmed the occurrence of su- 
percritical disk accretion onto BHs (Ohsuga et al. 2005, 
hereafter Paper I). 

In the case of the supercritical disk accretion flows 
around BHs, a large amount of photons generated inside 
the disk are swallowed by the BH with accreting matter. 
This effect largely attenuates the radiation force, assisting 
the occurrence of inflow motion. However, the situation 
considerably differs in disks around neutron stars (NSs). 
In this case, the radiation energy is not swallowed by the 
NS, and the energy of accreting matter would be converted 



to radiation energy at the NS surface. Thus, the radiation 
force is enhanced and works to prevent the inflow motion, 
although the matter is accelerated inward by the radiation 
force in cooperation with the gravity at the very vicinity 
of the BH. 

The X-ray source SS433 exhibits powerful relativistic 
jets. The observations of Doppler-shifted iron lines re- 
vealed that the mass-outflow rate of the jets highly ex- 
ceeds the critical rate for the NS/stellar BH (Kotani et al. 
1996). Thus, SS433 is thought to undergo supercritical 
accretion flows, although its central object is still not un- 
derstood. Whereas some researchers have suggested that 
a BH exists in SS433, a NS hypothesis is not denied. On 
the other hand, by the X-ray observations of low-mass X- 
ray binaries (LMXBs), Takahashi & Makishima (2006) re- 
cently discovered a sign of supercritical accretion onto the 
NSs. They reported that the spectra of LMXBs consist of 
three components in the supercritical phase. Two compo- 
nents are thought to be emission of the disk and the NS 
surface, which are observed when the X-ray luminosity is 
much smaller than the Eddington luminosity as well. The 
third component appears only in the supercritical phase, 
and is thus presumably the emission of outflows driven 
by the strong radiation force of the supercritical accretion 
disks. 

The super or near-critical accretion onto NS has been 
studied in spherical geometry by analytic method and nu- 
merical calculations (Burger & Katz 1980, 1983; Klein 
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et al. 1980; Miller 1990). They reported the maximum 
accretion rate to be about critical late and the maximum 
luminosity to be near the Eddington luminosity. However, 
Dcmmel et al. (1990) proposed a inhomogeneous spherical 
accretion model, by which the luminosity can exceed the 
Eddington luminosity for large mass-accretion rate. The 
maximum luminosity of supercritical disks around a mag- 
netized NS was suggested to be 46Le by Lipunov (1982). 
Two-dimensional RHD simulations of supercritical accre- 
tion flows around NSs were attempted by Okuda & Fujita 
(2000). However, though the matter was injected at the 
supercritical rate from the outer boundary, most of the 
matter was blown off from the system as the jets and the 
disk wind. The steady structure and physical properties 
of supercritical disk accretion onto a NS is still open issue. 

In this paper, we consider the supercritical disk accre- 
tion flows around non-magnetic and non-rotating NSs, 
based on long-term two-dimensional RHD simulations. 
We also compare them with supercritical BH accretion 
flows. In §2, our model and numerical method are de- 
scribed. We present the numerical results in §3. §4 and 
§5 are devoted to discussion and conclusions. 

2. MODEL AND METHOD 

2.1. Basic Equations 

We describe the RHD equations including the viscosity 
term. We use spherical polar coordinates (r,9,ip), where r 
is the radial distance, 6 is the polar angle, and ip is the az- 
imuthal angle. The origin is set at the central object (NS 
otherwise BH). In the present study, we assume axisym- 
metry with respect to the rotation axis (i.e., 8/ dip = 0), 
as well as reflection symmetry relative to the equatorial 
plane (with 9 = 7r/2). We describe the gravitational field 
of the central object in terms of pseudo-Newtonian hydro- 
dynamics, in which the gravitational potential is given by 
—GMJ (r — rg) as was introduced by Paczynsky & Wiita 
(1980), where M is the mass of the central object and 
rg = 2GM/c 2 is the Schwarzschild radius. The flow is as- 
sumed to be non-self-gravitating. The basic equations are 
the continuity equation, 
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the energy equation of the gas 
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and the energy equation of the radiation, 
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Here, p is the mass density, v — (v r ,vg,v v ) is the velocity, 
p is the gas pressure, e is the internal energy density of 
the gas, B is the blackbody intensity, Eq is the radiation 
energy density, F is the radiative flux, P is the radiation 
pressure tensor, k is the absorption opacity, q= (<lr,qg,q v ) 
is the viscous force, and <& v is is the viscous dissipative 
function, respectively. The radiation force / rad = (f r ,fe) 
is given by 
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where x(= K + P a T / Tn p ) is the total opacity with or being 
the Thomson-scattering cross-section and m p being the 
proton mass. 

As the equation of state, we use 

p=( 7 -l)e, (8) 

where 7 is the specific heat ratio. The temperature of the 
gas, T, can then be calculated from 

pkT 



P = 



fJ>m p 



(9) 



where k is the Boltzmann constant and p, is the mean 
molecular weight. 

The radiative flux and the radiation stress tensor are 
solved under the flux-limited diffusion approximation, so 
that they are expressed in terms of the radiation energy 
density (Levermore & Pomraning 1981). It gives correct 
relations in the optically thick diffusion limit and optically 
thin streaming limit, respectively. 

We assume that the only r</5-component of the viscous 
stress tensor, which plays important roles for the transport 
of the angular momentum and heating of the disk plasma, 
is non zero, 
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where 77 is the dynamical viscosity coefficient. Then, the 
radial and polar components of the viscous force are null 
(Sr = Qe =0), and the right hand side of equation (4) is 
described as 
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The viscous dissipative function is given by 
i 2 
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dr 



The dynamical viscosity coefficient is prescribed as 
P + A£ 



77 = a- 
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where a is the viscosity parameter, and Qk is the 
Keplerian angular speed. Here, the parameter, A, is given 
by (2 + TZ)/(6 + 3TZ + 1Z 2 ), using the dimensionless quan- 
tity, 1Z = \S7E \/ (x-Eq). Since it becomes 1/3 in the 
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optically thick regime and vanishes in the optically thin 
regime, the viscous stress tensor is related to the gas pres- 
sure or the total pressure, ptotai, as, 

f "Ptotai (optically thick limit) , . 

rtp 1 ap (optically thin limit) ' ^ ' 

It is modified a prescription of the viscosity, which was 
proposed by Shakura & Sunyaev (1973). 

2.2. Numerical Method 

We numerically solved the set of RHD equations us- 
ing an explicit-implicit finite-difference scheme on the 
Eulerian grids. Our methods were similar to those of 
Paper I, but we used the explicit scheme for solving the 
radiative flux term in equation (6), and employ constant 
time step (Ohsuga 2006). The computational domain 
consisted of spherical shells of 2.4rs < r < 500rs and 
<0 < 0.5-7T, which were divided into 96 x 96 grid cells. 

We employed two types of inner-boundary conditions 
for the accretion flows around the NSs (NS model) and the 
BHs (BH model) . The radial component of the radiative 
flux is given by 

F£ = -cE , (15) 

at the inner boundary for the BH model. This means that 
the radiative flux term always contribute to transport of 
the radiation energy onto the BH (Fq < 0). In addition, 
because of the inflow motion at the inner boundary (v r < 
0), the radiation energy and the gas energy are swallowed 
by the BH via advection, v r Eo < and (e + pv 2 /2)v r < 
0, in our simulations. In other words, the BH always 
swallows both the matter and the energy. In contrast, we 
assumed that the NS cannot swallow the energy. The gas 
energy advected onto the NS surface is supposed to be 
immediately released as the radiation. Hence, we set 

Fq = -ev r - ^pv 2 v r , (16) 

at the inner boundary for the NS model. Our simulations 
gave negative radial velocity at the inner boundary (v r < 
0), leading an outward radiative flux (Fq > 0) at the NS 
surface. When we updated the radiation energy density 
by solving the equation (6), the advective term of the 
radiation energy was numerically set to be null at the 
inner boundary, v t Eq = 0, whereas the radial velocity and 
the radiation energy density arc both non zero there. To 
sum up, the matter accretes onto the NS surface, but the 
energies of radiation and gas were not swallowed by the 
NS in our simulations. 

We employed a no-torque condition at the surface of the 
BH. Thus, the viscosity does not work at the inner bound- 
ary in the case of the BH accretion flows. In contrast, we 
assumed that the viscosity works even at the very vicin- 
ity of the inner boundary in the NS model. Therefore, 
the rotation energy of the gas is effectively dissipated and 
converted to radiation energy near the NS surface. This 
viscous effect at the inner boundary would depend on the 
rotation velocity of the NS. However, we note that our 
results do not change so much unless the rotation velocity 
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of the NS surface is comparable to the Keplerian velocity. 

The radiation energy density increases near the NS sur- 
face because of such inner-boundary conditions, leading 
the enhanced radiation force. The radiation force coun- 
teracts with the gravity and works to prevent the inflow 
motion. In contrast, the matter is accelerated inward by 
the radiation force in cooperation with the gravity in the 
case of the BH model. 

We started the calculations with a hot, rarefied, and op- 
tically thin atmosphere. There was no cold dense disk in 
the computational domain, initially. Matter was contin- 
uously injected into the computational domain through 
the outer disk boundary (r = 500rs, 0.45-7T < 6 < 0.57r). 
Therefore, we could avoid any influence of the initial con- 
figuration on the final result, although a long integra- 
tion time is required. The injected matter is supposed 
to have a specific angular momentum corresponding to 
the Keplerian angular momentum at r = 100rg. We set 
the mass- input rate, Mi nput , which is the mass injected 
into the computational domain per unit time, so as to be 
constant at the outer disk boundary. At the outer bound- 
ary region above the outer disk boundary (r = 500rs, 
9 < 0.457r), we used free boundary conditions and allowed 
for matter to go out but not to come in. 

Throughout the present study, we assume a = 0.1, 7 = 
5/3, (i = 0.5, and M = 1.4M . Thus, the radius of the 
inner boundary, 2.4rg, is nearly 10km. Such mass and the 
radius are consistent with recent models of a NS, although 
they slightly depends on the equation of state inside the 
NS. Since we aim to understand any differences of the 
accretion flows around NSs and BHs, the mass of the BH 
was also set to be 1.4M0. For the absorption opacity, we 
considered free-free absorption and bound-free absorption 
for solar metallicity (Hayashi, Hoshi, & Sugimoto 1962; 
Rybicki & Lightman 1979) 

3. RESULTS 

3.1. Quasi- steady Structure 

We first give a brief overview of the time evolution of 
our simulations. The overall evolution of the NS accre- 
tion flows is similar to that of the BH accretion flows 
(see Paper I in detail). It is divided into two distinct 
phases: an accumulation phase and a quasi-steady phase. 
Because the injected matter has a small specific angular 
momentum, it firstly free falls, but accumulates around 
regions with a radius of r ~ 100rs, forming a rotating disk 
(accumulation phase). Eventually, the viscosity starts to 
work. The angular momentum of the gas is transported 
outward, driving the inflow motion (quasi-steady phase). 
Radiatively driven outflow is generated around the rota- 
tion axis in this phase. The critical time separating these 
two phases roughly coincides with the viscous timescale. 

Figure 1 indicates time-averaged contours of the density 
in the quasi-steady state for the NS model (the upper-left 
panel) and the BH model (the upper-right panel), where 
we adopt Mi nput = 10 3 Le/c 2 . As shown in these pan- 
els, a high-density region appears around the equatorial 
plane (yellow) and a rarefied region forms near the rota- 
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Fig. 1. Time-averaged distributions of the density (the upper-left panel) and the radial velocity normalized by the speed of light 
(the lower-left panel) on the R-z plane for the NS model. Those for the BH model are shown in the upper-right and lower-right 
panels. The adopted parameter is M; npu t = 10 3 Le/c 2 . 



tion axis (blue) in both models. The high-density regions 
indicate the radiation pressure-supported disks. The ra- 
diation force drives the outflow in the less-dense region 
above the disk. We find a significant difference in the very 
vicinity of the inner boundary. Although the rarefied re- 
gion connects to the inner boundary along the rotation 
axis in the BH model, the high-density region (orange) 
surrounds the NS. 

Such a shell-shaped structure is built up by the radia- 
tion force, which counteracts with the gravity. Since the 
radiation energy, as well as the gas energy, is not swal- 
lowed by the NS in our simulations, the radiation energy 
increases near the inner boundary. Hence, the radial com- 
ponent of the radiative flux becomes positive (outward 
flux) there, so that the radiation force works to prevent 
the inflow motion (discussed later). Here, we note that 
such a difference in the density distribution would affect 



the emergent spectra and its viewing-angle dependence. 
We discuss this point in §3.4. 

In the lower-left (NS model) and lower-right (BH 
model) panels, we display the radial velocity normalized 
by the speed of light. It is found that the gas slowly falls 
inward in the dense disk regions. On the other hand, the 
matter is accelerated outward in the rarefied region above 
the disk, except for the very vicinity of the BH, where 
the gas falls onto the BH. We find in the NS model that 
the typical velocity of the outflow is 0.2 — 0.3c and its 
collimation angle is ~ 20°. The BH model exhibits more 
collimated high velocity regions along the rotation axis. 
The supercritical accretion onto the NS forms more pow- 
erful outflows than does that around the BH. Both the 
mass-outflow rate and the kinetic energy output rate by 
the outflow are larger in the NS model than in the BH 
model (discussed later). 
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Fig. 2. Time-averaged one-dimensional profiles of the den- 
sity (the top panel), the radial and rotation velocities (the 
middle panel), and the radiation temperature (the bottom 
panel) near the equatorial plane, 8 = 0.49-7T. The solid and 
dotted lines indicate the results of the NS model and the 
BH model, respectively. The Keplerian velocity is repre- 
sented by the dashed line in the middle panel, which almost 
overlaps with upper dotted line. The adopted parameter is 
Mmput = 10 3 i E /c 2 . 

In Figure 2 we represent one-dimensional distributions 
of the gas density (the top panel), the inflow and rota- 
tion velocities (the middle panel) , and the radiation tem- 
perature (the bottom panel) near the equatorial plane, 
9 = 0.497T. Here, the solid and dotted lines indicate the 
results of the NS model and the BH model, respectively. 
The radiation temperature is defined as T r = (Eq/o) 1 ^ 4 , 
where a is the radiation constant. As shown in the top 
panel, the gas density in the NS model is comparable to, 
or slightly smaller than, that in the BH model. It is found 
in the middle and bottom panels that there are significant 
differences in the profiles of the velocities and the radia- 
tion temperature between two models. These are induced 
by the differences in the inner-boundary conditions. 

In the present simulations, the viscosity worked even at 
the inner boundary for the NS model, whereas we used 
no-torque condition at the surface of the BH. Therefore, 
the rotation velocity deviates from the Keplerian velocity 
as the flow approaches the NS, though the gas rotates at 
about the Keplerian velocity in the regions of r>5rs (see 
the middle panel). In contrast, the rotation velocity of 
the BH model roughly agrees with the Keplerian velocity 
even in the vicinity of the inner boundary. 

We find in the bottom panel that the radiation energy 
is larger in the NS model than in the BH model near the 
inner boundary. Near the surface of the NS, the viscosity 
effectively converts the rotation energy into the internal 



energy of the gas. Subsequently, the gas-radiation inter- 
action turns such energy into the radiation energy In 
addition, the energy of the gas falling onto the NS is con- 
verted to radiation energy at the surface of the NS (see 
equation [16]). Due to the combination of these effects, 
the radiation energy increases around the inner boundary 
in the NS model. In the case of the BH accretion flow, the 
no-torque condition is employed and the energy is swal- 
lowed by the BH, attenuating the radiation energy density 
around the inner boundary. 

We find in the middle panel that the resulting inflow 
velocity — v r , in the NS model is much smaller than that 
in the BH accretion flows around the inner boundary. The 
radiation force is enhanced in the vicinity of the NS sur- 
face via the large radiation energy and its steep profile (see 
the bottom panel) . Thus, the radiation force works to de- 
celerate the inflow motion. That radiative deceleration in 
the supercritical flows around the NS was mentioned by 
Burger & Katz (1980; 1983). In the BH model, the radia- 
tive deceleration is relatively inefficient, since the flatter 
profile of the radiation energy density leads to a small 
radiative flux (see the bottom panel). In addition, the 
radiation force accelerates the gas inward in cooperation 
with the gravity at the surface of the BH, since the ra- 
dial component of the radiative flux is negative (inward 
flux). As a result, the inflow velocity increases as the flow 
approaches the BH. 

3.2. Mass Accretion Rate and Mass Outflow Rate 

We plot in Figure 3 the mass-accretion rate onto the 
NS/BH (circles) and the mass-outflow rate (squares) nor- 
malized by the critical accretion rate as functions of the 
mass-input rate. Here, the mass-outflow rate is the mass 
ejected through the (niter boundary per unit time by the 
high- velocity outflow whose radial velocity exceeds the es- 
cape velocity. The filled and open symbols indicate the 
results of the NS and BH models, respectively. 

We find that both the mass-accretion rate and the mass- 
outflow rate increase with an increase of the mass-input 
rate. We find that the supercritical accretion is not a 
peculiar phenomenon in BH accretion flows. The mass- 
accretion rate of the NS model highly exceeds the critical 
accretion rate, although it is as small as 20 — 30% for the 
BH model. Especially surprising is that the mass-outflow 
rate for the NS model exceeds that in the case of BH ac- 
cretion flows, although the mass-accretion rate is much 
smaller in the NS model than in the BH model. The su- 
percritical disk around NS blows away a few times as much 
matter as outflowing matter of the disk around the BH. 
The supercritical accretion flows around the NSs effec- 
tively produce the radiatively driven outflow. We discuss 
this point again in the next subsection. 

3.3. Luminosity and Kinetic Energy Output Rate 

In Figure 4 we represent the luminosity and the ki- 
netic energy output rate as functions of the mass-accretion 
rate onto the NS/BH, M acc . The filled (open) circles and 
squares indicate the luminosity, L ra d, and the kinetic en- 
ergy output rate, Lkini for NS (BH) model. Here, the 



Ohsuga 



[Vol. 



10 2 



10 s 



10 1 




m 0Ut (NS) 
m out (BH) 



m a =c(NS) 



10 s 



10 3 
M. ,c 2 /L v 

input / E 



io 4 



10 1 



~n r 




i- rad A E (BH) 



i kl „A E (BH) 



10 1 



10 £ 



10 E 



M acc c 2 /L E 



Fig. 3. Mass-accretion rate onto the NS/BH and mass-out- 
put rate by the high velocity outflows normalized by the crit- 
ical rate, m aC c (circles) and m ou t (squares), as functions of 
normalized mass-input rate. The filled and open symbols in- 
dicate the results of the NS model and the BH model, respec- 
tively. 

luminosity is evaluated by integrating the radiative flux 
at the outer boundary. The kinetic energy output rate 
is the kinetic energy ejected through the outer boundary 
per unit time by the high-velocity outflow. The luminos- 
ity of the NS, Lnsj which is evaluated by integrating the 
radiative flux at the inner boundary (see equation [16]), 
is represented by the triangles. The mass-input rate is set 
to be M input /(i E /c 2 ) = 3 x 10 2 , 10 3 , and 3 x 10 3 from left 
to right. 

We find in this figure that the luminosity and the ki- 
netic energy output rate increase with an increase of the 
mass-accretion rate in both models. It is also found that 
the energy conversion efficiency, (X ra d + / M acc c 2 , of 
the NS accretion flows is by one order of magnitude 
larger than that of the BH accretion flows. For in- 
stance, when the mass-accretion rate is around 70Le/c 2 , 
we find L ra d ~ 9Le and ikin ~ 19Ze for the NS model 
and L ra d ~ 2Le and Lkin ~ -^e for the BH model. The 
sum of the luminosity and the kinetic energy output rate 
almost equals the accretion energy for the NS model, 
L Ta d + -^kin ~ GMM acc /(2Ars — rg). In contrast, we find 
irad + ikm<GAfAf acc /(2.4r s -r s ) in the BH model. This 
is natural since, unlike a BH, a NS neither swallow the ra- 
diation energy nor the gas energy in the present study. 

In Figure 4 we find L ra d > L^ n for the BH model, but 
£rad < ^kin for the NS model. The radiation force is 
enhanced via the inner-boundary conditions for the NS 
model, whereby the rotation energy is converted to radi- 
ation energy via the viscosity and the energy of the in- 



Fig. 4. Luminosity (circles) and kinetic energy out- 
put rate (squares) as functions of normalized mass-ac- 
cretion rate onto the central object, M acc c 2 /LE, for 
Mn P ute 2 /LE = 3 X 10 2 , 10 3 ,3 X 10 3 from left to right. The 
filled and open symbols indicate the results of the NS model 
and the BH model, respectively. The luminosity of the NS is 
represented by the filled triangles. 

flowing matter is also converted to the radiation energy 
without being swallowed by the NS, driving the strong 
radiatively driven outflow. 

Next, we explain the reason why the mass- accretion rate 
onto the NS can be over the critical rate, L^/c 2 . In the 
present simulations, the luminosity of the NS corresponds 
to the energy of the accreting gas per unit time (see equa- 
tion [16]). Thus, if the gas sufficiently releases its en- 
ergy before reaching the inner boundary, the gravity is 
predominant over the radiation force at the NS surface, 
allowing the accreting motion. Otherwise, the radiation 
force would prevent the inflow motion at the inner bound- 
ary. Due to such mechanism, the resulting luminosity of 
the NS is sclf-rcgulated so as to meet the condition that 
the radiation force does not exceed the gravity. Since we 
employ the pseudo-Newtonian potential, the critical lu- 
minosity, in which the radiation force balances with the 
gravity, is ~ 3£e- To conclude, although the energy trans- 
ported onto the NS surface per unit time is limited to be 
^<3Le, the mass-accretion rate can greatly exceed the crit- 
ical rate, as long as the matter sufficiently loses its energy 
before reaching the NS surface. Since the BH can swal- 
low the gas energy, as well as the radiation energy, the 
mass- accretion rate can exceed the critical rate even if 
the matter has a large amount of the energy. 

3.4- Effective Temperature Profile 

Finally, we consider the effective temperature by solving 
radiation transfer, and represent its profiles for the inch- 
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nation angles of i = 0, 7r/12, and ir/6 in Figure 5. Here, 
X is the horizontal coordinate on the observer's screen. 
The adopted parameter is Mi npu t = 10 3 Le/c 2 . Then, the 
luminosity is ~ 5Le for the NS model and ~ 3Le for the 
BH model. In this subsection, we focus on the effective- 
temperature profile and its viewing-angle dependence. 

This figure shows that the maximal effective tempera- 
ture is not very sensitive to the viewing angle for the NS 
model. It merely varies from 3 x 10 6 K to 5 x 10 6 K. In 
contrast, it shifts from 2 x 10 6 K to 9 x 10 6 K for the BH 
model. We also find that the accretion flow around the 
NS exhibits a flatter profile of the effective temperature, 
even in the face-on view. 

As shown in the Figure 1, the density distribution is 
more anisotropic in the BH model than in the NS model 
(note that the outflow region above the disk is relatively 
dense in the NS model, as compared with the BH model). 
In addition, a high-density shell forms only in the NS 
model. Due to such enhanced anisotropy of the density 
distribution, the emission of the BH model is more colli- 
matcd in the polar direction than that of the NS model, 
producing the steeper effective-temperature profile in the 
face-on view (see also Paper I). In contrast, the effective- 
temperature profile becomes almost flat in the NS model. 
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Fig. 5. Effective-temperature profiles in the X-dircctions for 
the viewing angle of i = 0, 7r/12, and tt/6. Here, X is the 
horizontal coordinate on the observer's screen. The solid and 
dotted lines represent the results of NS model and BH model, 
respectively. 



4. DISCUSSION 

4-1. Comparison with SS433 jets 

By the present simulations, we showed that the rela- 
tivistic outflows appear above and below the supercritical 
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disk. Here, we compare the simulated outflows of the NS 
model with the jets of SS433. 

Based on the observations of the Dopplcr-shiftcd ion 
lines, Kotani et al. (1996) reported the mass-outflow 
rate of the SS433 jets to be > 5 x lO~ 7 M yr _1 , prob- 
ably 2 x 10 _6 M Q yr _1 . Since these values correspond 
to > 1.7 x 10 2 L E /c 2 and 6.7 x l0 2 L E /c 2 for the mass of 
1.4Mq, our NS model can explain the mass-outflow rate of 
the SS433 jets, if the mass-input rate is within the range 
3 x 10 2 <M inpu t/(L E /c 2 ) <3 x 10 3 . The typical velocity of 

the outflows, (2Lkin/moutiE) 1 ^ 2 c, is 0.2 — 0.3c. It nicely 
agrees with the velocity of 0.26c observed in SS433. 

However, the collimation angle of the jets of our model 
deviates from that of the SS433 jets. The collimation an- 
gle of the SS433 jets seems to be as small as a few degrees. 
However, we find the angle of ~ 20° in our simulations (see 
the upper-left panel in Fig. 1). The magnetic effects might 
resolve this issue. Kato, Hayashi, & Matsumoto (2004) 
have reported that the collimatcd jets are produced via 
a magnetic interaction between a weakly magnetized NS 
and its accretion disk. The jets might be collimated at 
the outer region, although the size of the computational 
domain was restricted to 500rs in the present study. We 
will clarify this point in future work by performing the 
simulations with larger computational domains. 

We next discuss about the density and the temperature. 
Kotani et al. (1996) also estimated the density and the 
temperature in the X-ray emitting region, ?^>10 9 cm (2 x 

10 3 r s for the mass of 1.4M ), as p < 10~ 8 gcm~ 3 and T ~ 
3 x 10 8 K. In Figure 6 we represent the resulting profiles 
of the density (the top panel) and the temperature (the 
bottom panel) near the rotation axis, 8 = 4° for Aii nput = 
10 3 i E /c 2 . As shown in this figure, the density and the 
temperature are roughly proportional to r^ 1 and r~ 2 / 3 , 
respectively. Thus, our results are roughly consistent with 
the observations as long as the slope of the profiles does 
not drastically change up to the regions of r ~ 2 x 10 3 rs, 
although our computational domain does not cover the X- 
ray emitting region. We need the numerical simulations 
with larger computational domains. 

The density and the temperature of the SS433 jets near 
the compact object are investigated by Inoue et al. (2001). 
They analytically studied and reported p^7x 10 _7 gcm~ 3 
at r~9x 10 7 cm. Since this size corresponds to r ~ 
2 x 10 2 rg for a mass of 1.4M , their result roughly agrees 
with our result (see the top panel). On the other hand, 
whereas we used the one-temperature assumption in the 
present simulations, it is not guaranteed in the outflow 
regions. This is because the equipartition time, <10 -3 
s, is comparable to, or larger than, the escape time at 
the regions of r^<200rs. In the case in which the one- 
temperature assumption breaks down, the ion tempera- 
ture would be much higher than the electron tempera- 
ture. Such a tendency was also reported by Inoue et al. 
(2001). In their work, the ion temperature was one or- 
der of magnitude higher than the electron temperature. 
However, their model predicts a much higher electron 
temperature, ~2x 10 10 K, than the result of our model, 
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Fig. 6. Time-averaged one-dimensional profiles of the den- 
sity (the top panel) and the gas temperature (the bottom 
panel) in the outflow regions (0 = 4°) for the NS model. 

~ 10 9 K. Although we solved a set of RHD equations in 
the present simulations, Inoue et al. (2001) investigated 
the only the energetics by assuming the dynamics and the 
shape of the jets. We should perform RHD simulations, 
in which the ion and electron temperatures are separately 
treated. 

4-2. Future work 

We need to calculate the spectra and compare them 
with the results of the observations for detailed study. 
A symptom of the supercritical accretion onto NSs has 
been discovered by the X-ray observations. Takahashi & 
Makishima (2006) reported that the LMXB spectra can 
be understood by the superposition of emissions from an 
accretion disk, the NS surface, and outflows. Since the 
emission, which is thought to come from the outflows, 
is not detected in the subcritical phase, they concluded 
that the outflows are accelerated by the strong radia- 
tion force of the supercritical accretion flows. By solving 
the monochromatic radiation-transfer equation, we calcu- 
lated the effective-temperature profiles, and show them 
in Figure 5. However, the emergent spectra of the su- 
percritical accretion flows are not a simple superposition 
of blackbody spectra with various effective temperatures 
(Ohsuga, Mineshige, & Watarai 2003). In addition, the 
inverse Compton scattering would also contribute to the 
spectra, though it is not included in our simulations. A 
detailed study of the emergent spectra should be explored 
in future work. 

In the present study, we consider the accretion flows 
around non-rotating NSs. Our results do not change that 
much as long as the rotation velocity of the NS surface 
is much smaller than the Keplerian velocity. Here, we 
note that the spin-up timescale of the NSs by mass accre- 
tion, ~ 10 8 yr, is much longer than the integration time 
of the present simulation, ~ 10 s. The rapidly rotating 
NSs seem to exist in the LMXBs, since they are thought 
to be progenitors of millisecond radio pulsars (Backer et 
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al. 1982; Wijnands & van der Klis 1998). The rapid ro- 
tation plays important roles in the disk structure and the 
temperature profile via modification of the space-time ge- 
ometry (Bhattacharyya et al. 2000, 2001a, b). A detailed 
consideration of the disks around the rapidly rotating NSs, 
however, beyond the scope of this paper. 

We consider the non-magnetized NSs in the present 
study. If the NSs have strong magnetic fields, the flow 
structure might drastically change (Kato et al. 2001; 
Koldoba et al. 2002; Romanova et al. 2003). In addi- 
tion, the source of the disk viscosity is likely to be of 
magnetic origin, though we used the a-prescription of 
the viscosity (Hawley, Balbus, & Stone 2001; Machida, 
Matsumoto, & Mineshige 2001; Balbus 2003, for a re- 
view). We need radiation magnetohydrodynamic (MHD) 
simulations. Local radiation MHD simulations of the ac- 
cretion flow have been performed by Turner et al. (2003) 
and Turner (2004). Furthermore, we should include the 
synchrotron sclf-Compton processes. It is important to 
separately treat the ion temperature and the electron tem- 
perature. As we discussed in §4.1, the one-temperature as- 
sumption might break down in the outflow regions because 
of low density and high velocity. The radiation MHD sim- 
ulations are required again for such a study. 

Although the matter did not lose angular momentum by 
the gas-radiation interaction in the present study, the ra- 
diation drag would extract angular momentum from the 
disk. However, we here stress that our results do not 
change very much even if we take the radiation drag into 
consideration. Miller & Lamb (1993) have reported that 
the mass accretion is enhanced via radiation drag by the 
radiation from the NSs. The radiation drag- induced mass- 
accretion rate is given by Lns/c 2 , if all photons emit- 
ted at the NS surface contribute to extract the angu- 
lar momentum (sec also Umcmura 2001). This value is 
much smaller than the resulting mass-accretion rate (see 
Figure 4) . Additionally, the radiation drag- induced mass- 
accretion rate by the radiation from the disk is L ra d/c 2 
maximum. It is also smaller than our results as shown in 
Figure 4. Since photons are emitted from rapidly rotating 
matter, and since some part of the photons escape through 
the rarefied region above the disk without any interaction 
with the disk matter, the efficiency of the radiation drag 
would decrease. Thus, we conclude that our present re- 
sults do not change very much due to radiation drag. For 
a detailed study, we should sclf-consistently solve the ra- 
diation fields and the accretion flow, by performing the 
rclativistic RHD simulations. 

5. CONCLUSIONS 

By performing two-dimensional RHD simulations, we 
investigated the supercritical accretion flows around the 
NSs. We summarize our results as follows. 

(1) Since the radiation energy, as well as the gas energy, 
is not swallowed by the NS, and since the rotation energy 
is effectively converted to radiation energy at the inner- 
most regions via the viscosity, the radiation force is en- 
hanced near the NS surface. The enhanced radiation force 
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counteracts with the gravity and produces high-density 
regions around the NSs. Such a high-density shell is a 
peculiar structure of supercritical accretion flows onto the 
NSs, although the disk and outflows are also produced in 
the case of the BH accretion flows. 

(2) Although the enhanced radiation force reduces the 
mass-accretion rate, and drives the strong outflows above 
and below the disk, supercritical disk accretion onto the 
NS is possible. The mass-accretion rate onto the NS is 
20 — 30% of that onto the BH, on the condition that we 
employ the same mass-input rate. The mass-outflow rate 
is a few-times larger in the NS accretion flows than in the 
BH accretion flows. 

(3) The accreting matter releases most of its energy be- 
fore reaching the NS surface. Due to this mechanism, the 
luminosity of the NS is insensitive to the mass-accretion 
rate. Therefore, although the matter accretes onto the 
NS at the supercritical rate, the radiation force does not 
exceed the gravity at the NS surface. 

(4) The supercritical accretion flows onto the NSs 
mainly release the accretion energy as the kinetic energy 
of the outflows, especially for a large mass-accretion rate. 
In contrast, the kinetic energy output rate via the out- 
flows does not exceed the disk luminosity in the case of 
BH accretion flows. 

(5) Our simulations succeed in reproducing the observed 
properties of the SS433 jets, except for the collimation 
angle. The resulting jet velocity (0.2 — 0.3c) and mass- 
outflow rate (150 — 700Le/c 2 ) agree with the observations. 
However, the collimation angle of the jet in our simula- 
tions (~ 20°) is larger than that of SS433 (a few degrees). 
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